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As they decelerate through the atmosphere, meteors deposit mass, momentum and en- 
ergy into the surrounding air at tremendous rates. Trauma from the entry of such bolides 
produces strong blast waves that can propagate hundreds of kilometers and cause substan- 
tial terrestrial damage even when no ground impact occurs. We present a new simulation 
technique for airburst blast prediction using a fully-conservative, Cartesian mesh, finite- 
volume solver and investigate the ability of this method to model far-field propagation 
over hundreds of kilometers. The work develops mathematical models for the deposition of 
mass, momentum and energy into the atmosphere and presents verification and validation 
through canonical problems and the comparison of surface overpressures, and blast arrival 
times with actual results in the literature for known bolides. The discussion also examines 
the effects of various approximations to the physics of bolide entry that can substantially 
decrease the computational expense of these simulations. We present parametric studies 
to quantify the influence of entry-angle, burst-height and other parameters on the ground 
footprint of the airburst, and these values are related to predictions from analytic and 
handbook-methods. 
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I. Introduction 


N mid February 2013 an asteroid measuring approximately 20 meters in diameter entered the sky over 
Chelyabinsk Russia. The bolide had a mass of about 12,500 metric tons and vaporized nearly completely 
above 25 km. Despite its relatively small size, the asteroid carried with it tremendous kinetic energy and 
deposited approximately 520 kilotons into the atmosphere as it entered and burned up.! Over a minute 
later, the shock wave from this sudden deposition of energy reached the ground breaking glass and damaging 
structures in a region covering more than 5000 square kilometers. 
Government infrasound monitoring records an average of around 27 encounters with objects larger than 
a meter in diameter annually.2 Each decade, approximately seven encounters occur which release over 
10 kilotons into the atmosphere.*? Since the population of potentially hazardous objects follows a rough 
power law, the best current estimates predict encounters with Chelyabinsk-sized objects approximately every 
80 years.4 
In response to this threat, NASA’s Near-Earth Object program initiated a new research activity actively 
focused on Planetary Defense in October of 2014.° This program is structured around four thrusts: (1) Char- 
acterization and Composition, (2) Entry/Break-up physics, (3) Atmospheric Propagation and Impact Effects, 
and (4) Physics-Based Risk Assessment. 
The simulation and modeling effects discussed in the proposed paper target the third of these tasks, 
namely propagation of the airburst through the atmosphere and estimation of the impact effects. Clearly 
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this element is closely connected to both entry/breakup and risk assessment. In a nutshell, atmospheric 
propagation takes inputs from the entry modeling task and computes the land/water footprint to feed 
risk assessment. From a simulation and modeling standpoint, modeling entry physics typically involves 
highly-detaled multi-physics simulations done in the frame of the asteroid as it enters,° while atmospheric 
propagation is usually done in the Eulerian frame of a fixed observer on the ground and taking energy- 
deposition as input. 

In the literature, surface impact affects are often assessed using either analytic models or by model- 
ing the energy release through correlation with static nuclear detonations.“? More detailed analysis has 
also been performed using numerical simulation of the atmospheric blast propagation using various hy- 
drocodes.'3:® 19.11 Such simulations are often driven by either triggered breakup of a known object or 
by an estimate of the energy deposition profile which can be generated for specific bolide entries through 
observation of the light-curve,! or infrasound measurements.* !? 

The approach taken in the proposed paper uses the energy deposition profile as input and performs blast 
propagation using a fully-conservative finite-volume method with adaptive meshing. Rather than reconstruct 
a specific entry event, our goal in this work is to drive physics-based risk assessment. To this end, we desire the 
flexibility to drive simulation from outputs of entry-modeling task, analytic or PDE-based models,”® +!%"!° 
observed light curve data, or infrasound measurements. 


The final paper will include a theoretical section showing how energy deposition profiles are used to 
determine initial conditions for simulations taking energy deposition profiles as input. Details of the time- 
dependent modeling approach will be presented along with a discussion of computational expense. It will also 
include verification and validation data using a spherical blast model at relevant conditions. Investigations 
wil include examination of the influence of entry-angle, burst-height and other parameters. Further valida- 
tion will include comparisons with data from the Chelyabinsk meteor examining both ground over-pressure 
and blast arrival time using data from ground reports and tume-synched video. Notional titles for sections of 
the final paper are shown below along with some preluminary results for some of the investigations. 


II. Theoretical Development and Numerical Method 


A. Conservation of Mass, Momentum, and Energy 


During its violent entry, a meteoroid transfers mass, a 
momentum and energy into the surrounding air at 
tremendous rates. Observational data and models a 

of the physical processes are frequently processed to x 


express the net effects of the near-field physics as en- 
ergy deposition profiles in terms of energy-released 
as functions of both time and altitude. From the 
point of view of a fixed computation in an Eulerian 
frame on the ground, this deposition can be mod- 
eled as a time-dependent source for the state of the a 
fluid in a control volume surrounding the trajectory .~ " 
of the entering body. 
To aid in analysis, Figure 1 shows a tube-like Figure 1. Sketch showing control volume around tra- 
control volume sketched around the trajectory of Jectory of entering body used for converting energy- 


i h bh th h eh deposition profile data into time-dependent source 
i pony ero : oe the atmosphere at 1g terms for fluid state used in numerical simulations. 
speed.... 


a 


... Include theoretical development to support derivation of tume-dependent source terms for state variables 
used to drive simulations from energy-deposition profiles... 


B. Numerical Method 


Simulations for this work use a solver based on the Cartesian cut-cell approach of Reference [16] in which the 
governing equations are discretized on a multilevel Cartesian mesh with embedded boundaries. The mesh 
consists of regular Cartesian hexahedra everywhere, except for a layer of body-intersecting cut-cells which 
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intersect the geometry at boundaries. The spatial discretization uses a second-order accurate finite volume 
method with a weak imposition of boundary conditions, resulting in a system of equations 


R(Qr) =0 (1) 


The flux-vector splitting approach of van Leer!’ is used for residual evaluation. Although it consists 
of nested Cartesian cells, the mesh is viewed as an unstructured collection of control volumes making the 
approach well-suited for solution-adaptive mesh refinement. Steady-state flow solutions are obtained using 
a five-stage Runge-Kutta scheme with local time stepping and multigrid. Domain decomposition via space- 
filling curves permits parallel computation; for more details see Aftosmis et al. and Berger et al.1®"?° 


.. report modification of basic flow solver to include buoyancy effects through the Froude number and 
discussion of atmosphere model used in numerical investigations... 


C. Meshing and Computational Resources 


III. Numerical Investigations 
A. Verification for a Strong Spherical Blast 


... Report verification results for numerical investigation of spherical blast using a spherical charge of known 
strength compared with analytic solution of a equivalent blast emanating from a point source. Include a 
comparison of predicted overpressure with results from Friedlander equation. Demonstration to show numer- 
wcal method predicts well-known logarithmic asymptotic deceleration of blast wave as the radius of the blast 
expands... 


B. Chelyabinsk Meteor 


...Present validation of sumulation methodology using data collected for Chelyabinsk Meteor. Use energy 
deposition profile from Brown et al.? and compare with data for ground overpressure and shock arrival-time 
found in Reference [1]... 


Figure 2 shows simulation results for Chelyabinsk meteor using deposition profile from Reference [3]. 
These simulations used an entry angle of 18°, with peak brightness near the half-energy height at 29.5 km 
altitude and a total energy of 520 kt. The computational domain was 80 km high and extended 256 km in 
both directions along the ground plane. The simulation included a bit over 100 km of the terminal flight 
path. The computational mesh used ~90 M cells with roughly 1 m resolution along the trajectory and ~10 m 
on the ground in the region of peak overpressure. 

Figure 3 shows the ground overpressure footprint predicted by this simulation with an overlay of the 
entry trajectory as seen from above. The origin of the plot on the left is at the location of peak brightness. 
Isobars in both frames are have an increment of 1% of ambient pressure (sea level, standard atmosphere). 
The map at the right shows the predicted ground-footprint overlaid onto glass damage data reported by 
Popova and Jenniskens et al.' Statistical correlations in Mannan and Lees?! for blast-wave data indicate 
that an overpressure of ~700 Pa (0.69%) is able to shatter about 5% of typical glass windows while Ap of 
~6% breaks roughly 90%. The predicted overpressure contours shown in Figure 3 show very good correlation 
with glass damage reports. 

Figure 4 shows predicted shock-wave arrival time at ground level for comparison to data in Reference [1]. 
The time increment between successive contours is At = 33.2 seconds. Peak brightness is at the origin. The 
map in the frame at the right shows the ground track for correlating with other simulations and collected 
data. 


C. 90° Entry Angle 
Risk assessment is generally performed using monte-carlo simulations of static airbursts varying both the 


altitude and strength of the detonation.?? This allows handbook methods, like those in Reference [9], to be 
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Figure 2. Simulation results for Chelyabinsk meteor using deposition profile from Brown et al.? 18° entry 
angle, with peak brightness at 29.5 km altitude. Total Energy = 520 kt, 
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Figure 3. Ground footprint for Chelyabinsk meteor showing computed peak over-pressures overlaid on ground- 
damage map of Reference [1] using the energy deposition profile from Brown et al.* Increment between isobars 
is 1% of sea level pressure for standard atmosphere. Peak brightness is at (0,0). Statistical correlations in 
Mannan and Lees?! for blast-wave data indicate that an overpressure of ~700 Pa (0.69%) is able to shatter 
about 5% of typical glass windows while Ap of ~6% breaks roughly 90%. 
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Figure 4. Ground footprint showing computed shock-wave arrival time at ground level. Map at right shows 
ground-track. Peak brightness is at (0,0). 
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used to predict ground damage with statistical correlations. As we transition our risk assessment to physics- 
based blast predictions, the 90° entry case is of particular interest. Symmetry dictates that at this entry 
angle, the blast footprint will consist of a series of concentric circles describing variation in peak overpressure 
with distance from ground-zero. Since the energy deposition profile is roughly logarithmic with altitude, 
the vast majority of the energy along the trajectory is released near the point of peak brightness. Given 
these observations, its interesting to compare the ground overpressures from simulations along a 90° entry 
trajectory with those from a spherical detonation of the same magnitude at the same altitude. This case 
provides an important bridge between the detailed simulation of entry trajectories with more traditional 
static airburst correlations. 

Figure 5 presents a comparison of blast evolution from a line source with a Chelyabinsk-like deposition 
profile with that of a spherical charge with a 1 km radius for a 90° entry. The burst height for the spherical 
detonation was set to 29.5 km which corresponds closely to the half-energy altitude of the line source. Both 
simulations released a total of 520 kt to facilitate ground pressure comparisons with the 18° entry case 
presented earlier. Figure 5 presents the temporal evolution of both scenarios through several snapshots. 
The uppermost snapshots are taken shortly after energy release, at t* = 3.14. At the left, buoyancy-driven 
vortical structures travel upward along the trajectory while momentum in the core continues to advance 
toward the ground. At the right, the blast wave has propagated a radial distance of about 10 km from the 
blast center and the bubble has begun to flatten due to buoyancy effects acting upon the hot gas within. 
As the simulations advance in time, marked similarities begin to appear in the structure of the two blasts, 
and while the trajectory simulation retains a more intricate structure, its clear from both the color map 
and the location of the blast front that the strength and speed of the shock advancing toward the ground 
is remarkably similar in the two simulations. The final frames show the blast evolution at non-dimensional 
times of ~34.5 (roughly 114 seconds after detonation). Both blasts have reached the ground plane, and the 
eround footprint has expanded to a radius of ~30 km from ground zero in both simulations. 

Figure 6 shows a detailed comparison of the overpressure-envelope on the ground. For reference, the 
frame at the left shows a snapshot through the domain (via overpressure contours, (p — Poo) /Poo) taken at 
around 139 seconds after energy release. The right frame shows the envelope of recorded over /underpressure 
as well as instantaneous pressure profiles for both scenarios at 139 seconds. The yellow-shaded region shows 
the envelope for the spherical charge while the red line shows that of the line-source. Away from ground-zero, 
the envelopes for the two cases are remarkably similar, in particular the decay of the envelope with increasing 
radius is very well predicted by the spherical approximation. Nearer to the origin, the spherical charge is over 
predicts the line source however this appears to be an artifact of the particular energy deposition profile used 
to initiate this simulation which included a weak final fragmentation event near the end of the trajectory. 
Returning to Figure 5, evidence of this event is clearly visible in the top three frames of the blast from the 
line source and it appears to weaken the main blast front as it advances toward the ground near the origin. 
Since the spherical charge was initially homogeneous, it contains no such structures. 

The instantaneous pressure profiles shown by the red and blue dashed lines in the frame at the right of 
Figure 6 show remarkable similarities. At 139 seconds, both have expanded to a radius of just over 40 km. 
In addition to virtually identical propagation speeds, the two profiles show remarkably similar peak over- 
and under pressures. 


... Include follow-on discussion comparing results for spherical charge with overpressure and ground dam- 
age predictions using military-handbook correlations found in Reference [9]....... 


IV. Summary 


Me This work presented a new simulation technique for bolide-entry modeling using a fully-conservative, 
Cartesian mesh, finite-volume solver. The presentation included development of mathematical models for 
the deposition of mass, momentum and energy into the atmosphere and included studies for verification 
and validation of these models through canonical solutions of spherical blasts and the comparison of surface 
overpressures, and blast arrival time with actual results in the literature for the 2018 Chelyabinsk meteor. 

The discussion also examined various approximations to the physics of bolide entry to decrease the com- 
putational expense while retaining sufficient fidelity to drive ground damage estimates required for credible 
risk assessment. We present parametric studies that quantify the influence of entry-angle, burst-height and 
other parameters on the ground footprint of the airburst, and these values are related to predictions from 
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Figure 5. Comparison of blast evolution from a line source with a Chelyabinsk-like deposition profile with 
that of a spherical charge with a 1 km radius for a 90° entry. Burst height = 29.5 km, Total Energy = 520 kt, 
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Figure 6. Comparison of overpressure envelope from a line source with a Chelyabinsk-like deposition pro- 
file with that of a spherical charge with a 1 km radius for a 90° entry. Burst height = 29.5 km, Total 
Energy = 520 kt. Dashed lines show instantaneous ground pressure profiles shown for both simulations at 
139 sec after energy release. 


analytic and handbook-methods. .... 
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